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ABSTRACT 

We study the relation of AGN accretion, star formation rate (SFR), and stellar mass (M*) us¬ 
ing a sample of « 8600 star-forming galaxies up to z=2.5 selected with Herschel imaging in 
the GOODS and COSMOS fields. For each of them we derive SFR and M*, both corrected, 
when necessary, for emission from an active galactic nucleus (AGN), through the decom¬ 
position of their spectral energy distributions (SEDs). About 10 per cent of the sample are 
detected individually in Chandra observations of the fields. For the rest of the sample we 
stack the X-ray maps to get average X-ray properties. After subtracting the X-ray luminosity 
expected from star formation and correcting for nuclear obscuration, we derive the average 
AGN accretion rate for both detected sources and stacks, as a function of M*, SFR and red- 
shift. The average accretion rate correlates with SFR and with M*. The dependence on SFR 
becomes progressively more significant at z>0.8. This may suggest that SFR is the original 
driver of these correlations. We find that average AGN accretion and star formation increase in 
a similar fashion with offset from the star-forming “main-sequence”. Our interpretation is that 
accretion onto the central black hole and star formation broadly trace each other, irrespective 
of whether the galaxy is evolving steadily on the main-sequence or bursting. 
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1 INTRODUCTION 

A causal connection between super massive black hole (SMBH) 
and galaxy growth has been suggested hy a number of studies, 
based on empirical correlations between black hole mass and inte¬ 
grated galaxy properties: galaxy bulge M*, velocity dispersion (e.g. 
Magorrian et al. 1998; Gebhardt et al. 2000; Ferrarese 2002; Gul- 


* Herschel is an ESA space observatory with science instruments provided 
by European-led Principal Investigator consortia and with important partic¬ 
ipation from NASA, 
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tekin et al. 2009). In addition, the cosmic star formation history and 
the black hole accretion history follow parallel evolutionary paths, 
peaking at z~2 and declining towards the local Universe (Boyle & 
Terlevich 1998; Shankar et al. 2009; Madau & Dickinson 2014). 

Despite the mutual dependence on a common cold gas supply, 
such connections are not trivial given the vastly different spatial 
scales at which star formation (many kpc) and SMBH accretion 
(sub-pc) typically operate. Different scenarios have been proposed 
to justify the necessary loss of gaseous angular momentum, such as 
nuclear bars, minor and major merger events (e.g. Garcia-Burillo 
et al. 2005). However, the detailed mechanisms responsible for trig¬ 
gering black hole accretion and star formation are still poorly un- 
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derstood (e.g. see comprehensive review by Alexander & Hickox 

2012 ). 

Recent studies have highlighted a two-fold galaxy evolution¬ 
ary scheme. More than 95 per cent of star-forming galaxies fol¬ 
low a reasonably tight relation frequently called the “star-formation 
main-sequence” (MS, scatter is about 0.2-0.4 dex) between SFR 
and Mt, from the local Universe up to z~3 (Noeske et al. 2007; 
Daddi et al. 2007; Magdis et al. 2010; Elbaz et al. 2011; Whitaker 
et al. 2012; Schreiber et al. 2014; Speagle et al. 2014). This trend 
is currently thought to reflect a large duty cycle of steady star for¬ 
mation in galaxies, fueled by a continuous gas inflow (Dekel et al. 
2009; Ciotti et al. 2010). The most massive galaxies have larger 
gas reservoirs (Tacconi et al. 2013) and thus higher SFR. How¬ 
ever, there are a few (<5 per cent) outliers in the SFR-M* plane 
with > 4 times larger specific SFR (sSFR^). These off-sequence 
“starbursts” play a minor role in the cosmic star formation his¬ 
tory (Rodighiero et al. 2011) and show disturbed morphologies, 
probably due to galaxy interactions or gas-rich major mergers (e.g. 
Hopkins & Hernquist 2009; Veilleux et al. 2009). This distinction 
is supported by several studies, claiming a systematic variation of 
several galaxy properties with offset from the main-sequence; the 
off-sequence galaxies show more compact structures (Elbaz et al. 
2011; Wuyts et al. 2011b), warmer interstellar dust (Magnelli et al. 
2014), higher gas-to-M* ratio (Gao & Solomon 2004), larger Ear- 
to-Mid Infrared flux ratio (Nordon et al. 2010, 2012) and higher 
star formation efficiency (SFE^, Daddi et al. 2010b; Genzel et al. 
2010). The question of whether the transition from MS to starburst 
galaxies is steady or discontinuous remains open. 

A similar two fold-scheme is found also for the star-forming 
properties of low-luminosity X-ray-selected AGN (Lx <10"^^ erg 
s“^), showing at best a weak correlation between Lx and SFR, 
while bright Quasars follow a positive correlation with SFR at least 
up to z~L probably driven by major mergers (Lutz et al. 2008; 
Netzer 2009; Shao et al. 2010; Lutz et al. 2010; Rosario et al. 
2012). At z~2 such a correlation seems weak or absent (Rosario 
et al. 2012; Harrison et al. 2012). However, compelling evidences 
of AGN-driven feedback (e.g. Farrah et al. 2012) and the inverted 
correlation found in smaller samples of luminous z~2 AGN (Page 
et al. 2012) caution that our current picture of AGN/galaxy coevo¬ 
lution is still dependent on sample statistics and selection biases. 

By exploiting large samples of X-ray selected AGN, Mul- 
laney et al. (2012b) found that about 80 per cent of X-ray AGN 
live in main-sequence galaxies, 5-10 per cent in starburst galax¬ 
ies, and about 10-15 per cent in quiescent systems. While Santini 
et al. (2012) reported larger mean SFRs for the hosts of X-ray AGN 
compared to a mass-matched inactive reference that includes both 
star-forming and passive systems, Rosario et al. (2013) found very 
similar SFRs in X-ray AGN hosts and a mass-matched reference 
of only star-forming galaxies. All three studies reinforce the idea 
that most of the SMBH accretion is taking place in star-forming 
systems. 

In contrast to the studies of AGN hosts and to reach a compre¬ 
hensive understanding of the cosmic SMBH growth, several recent 
studies take a census of AGN accretion history on the basis of Far 
Infrared (FIR) and/or mass-selected samples of galaxies. Unlike the 
weak or absent correlation between star formation and black hole 


^ sSFR is defined as the ratio between SFR and M*. 

^ SFE is defined as the ratio between SFR and cold gas mass. 


accretion rate (BHAR^) for samples of AGN hosts, there is clear 
correlation of average BHAR and key properties of galaxy sam¬ 
ples. Positive and close to linear correlation is found between av¬ 
erage BHAR and mean stellar mass at various redshifts (Mullaney 
et al. 2012a) as well as with SFR (Rafferty et al. 2011; Chen et al. 
2013). These apparently contradictory results have been interpreted 
by Hickox et al. (2014) as due to different variability time scales be¬ 
tween nuclear activity and global star formation. According to this 
scenario, both components are intimately connected at any time: 
while star formation is relatively stable over ~100 Myr, the AGN 
might vary over ~5 orders of magnitude on very short (about 10® 
yr) time scales (Hickox et al. 2009; Aird et al. 2012; Bongiorno 
et al. 2012). In this scenario, all episodes of star formation are ac¬ 
companied by SMBH growth, but only when smoothing over the 
variations of individual sources do the average properties of AGN 
and their hosts show a consistent evolution, as stated by Mullaney 
et al. (2012a) and Chen et al. (2013). 

Both these latter studies derived average trends by binning 
their parent samples as a function of M* or SFR. Since their selec¬ 
tion techniques were mostly sensitive to main-sequence galaxies, 
in principle the resulting correlations found with average BHAR 
might be primarily due to one parameter, but reflected into a cor¬ 
relation with the other one, simply because of the main-sequence 
relation that holds between the two. To break this degeneracy and 
investigate in detail the role of AGN accretion in the context of 
galaxy evolution, it is necessary to split the sample as a function 
of both SFR and M* and study the evolution of the average AGN 
accretion properties in the SFR-M* plane at different redshifts. 

The primary goal of this work is to map the average BHAR 
as a function of SFR, M* and redshift. Our analysis exploits one of 
the widest compilations of FIR selected galaxies at non-local red¬ 
shifts. Robust SFRs for each individual source of the sample are 
measured from data taken by the Herschel Space Obsenatory (Pil- 
bratt et al. 2010). Our sample spans about three orders of magnitude 
in M* and four in SFR in the redshift range 0<z<2.5. For the first 
time we also investigate the role of AGN activity in off-sequence 
galaxies with respect to their main-sequence counterparts, seeking 
to constrain the parameter that primarily drives the growth of active 
SMBHs. 

We used the FIR data in COSMOS (Scoville et al. 2007) 
and from the Great Observatories Origins Deep Survey (GOODS) 
South (GOODS-S) and North (GOODS-N) fields, obtained with the 
//errcAeZ-Photodetector Array Camera and Spectrometer (PACS, 
Poglitsch et al. 2010), as part of the PACS Evolutionary Probe 
(PEP*', Lutz et al. 2011) project. In the GOODS fields, PEP data 
are also combined with the deepest observations of the GOODS- 
Herschel (GOODS-H®; Elbaz et al. 2011) open time key program. 
In addition, PACS observations at 70 (in GOODS-S only), 100 and 
160 pm are supplemented with sub-millimeter photometry at 250, 
350 and 500 pm obtained by the Spectral and Photometric Imag¬ 
ing Receiver (SPIRE, Griffin et al. 2010), as part of the Herschel 
Multi-tiered Extragalactic Survey (HerMES®, Oliver et al. 2012). 

The paper is structured as follows. In Section 2 we present our 
parent sample and multi-wavelength photometry. In Section 3 we 
introduce individual Herschel sources in the SFR-M* plane, while 


The tenns BFIAR and (Mnh) adopted throughout the paper are assumed 
to have the same physical meaning. 

* http : //www2011 . mpe . mpg . de/ir/Research/PEP/ 

® http : //hedam . oamp . fr/GOODS-Herschel 
® http://hermes.sussex.ac.uk 
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their average X-ray properties are derived in Section 4. We present 
the observed relationships between AGN and galaxy properties in 
Section 5, and discuss the implication of this work in Section 6. 
We list our concluding remarks in Section 7. Throughout this pa¬ 
per, we assume a Chabrier (2003) initial mass function (IMF) and 
a flat cosmology with 12^ = 0.30, Ga = 0.70, and Ho — 70 
kms“^ Mpc“^. 


2 SAMPLE DESCRIPTION 

Our sample exploits Herschel-PACS (Data Release 1) and SPIRE 
observations in the GOODS-North (GN hereafter), GOODS-South 
(GS hereafter) and COSMOS fields, covering in total about 2 deg^. 
In the following Sections we present the parent sample and briefly 
mention the cross-match with multi-wavelength identifications, re¬ 
ferring the reader to Lutz et al. (2011); Berta et al. (2011); Oliver 
et al. (2012) and Magnelli et al. (2013) for a detailed description of 
data reduction and construction of multi-wavelength catalogues. 

2.1 Ear-Infrared selected galaxies 

The parent sample includes all sources in the GOODS and COS¬ 
MOS fields wifh >3cr flux density in at least one PACS band. 
In both GOODS fields, FIR data are taken from the blind cata¬ 
logues described in Magnelli et al. (2013) that combine the data 
of PEP (Lutz et al. 2011) and GOODS-H (Elbaz et al. 2011). In 
the GOODS fields, flux density limits (3cr) reach « 1-2 mjy in the 
PACS bands, depending on filter and depth of observation, and « 8 
mJy in SPIRE-250 /rm. The confusion limit reachable with PACS 
ranges between 1.3 and 5 mJy (5 (j, Berta et al. 2011), while SPIRE- 
250 /rm observations are fully limited by confusion noise (Oliver 
et al. 2012). In the COSMOS field, the depth achieved by SPIRE- 
250 fim observations is more comparable to that reached by PACS 
ones, being « 5, 8 and 10 mJy at 100, 160 and 250 ^m, respec¬ 
tively. This is the main reason why we only used PACS-selected 
galaxies in the GOODS fields, while in COSMOS our selection also 
exploits SPIRE-selected ones. We note that the background level 
of PACS and SPIRE observations is relatively flat in each field. 
We checked fhat galaxies faken from the SPIRE-selected catalogue 
have consistent SFR and M* values (within a factor of two) with 
those taken from the PACS-selected sample at the same redshift. 
The presence of SPIRE-selected sources allows to double our sam¬ 
ple of star-forming galaxies without introducing a significant bias 
in our analysis. The extraction of PACS flux densities in all fields 
was performed blindly as described in Berta et al. (201 1), Lutz et al. 
(2011) and Magnelli et al. (2013), while for SPIRE sources it fol¬ 
lows the approach presented by Roseboom et al. (2012). 

2.2 Multi-wavelength identification 

The cross-match between PEP data and the extensive broad¬ 
band photometry available from the ultraviolet (UV) to the sub¬ 
millimeter has been accomplished via a maximum likelihood al¬ 
gorithm (Sutherland & Saunders 1992; Ciliegi et al. 2001) to 
deep Multiband-Imaging Photometer (MIPS) Spitzer detections at 
24 pm (Magnelli et al. 2011, 2013), whose positions have been 
used as priors to extract SPIRE fluxes in the sub-mm (Roseboom 
et al. 2012). The cross-match to optical/UV wavelengths in both 
GOODS fields is described in detail by Berta et al. (2010, 2011). 
We collected 892 (GS) and 850 (GN) FIR-selected sources with 
broad-band photometry from the optical/UV to the sub-mm. 


An extensive photometric coverage is available also in COS¬ 
MOS, where fluxes from the PEP catalogue have been cross- 
matched to 24 pm data (Le Floc’h et al. 2009), in turn used as po¬ 
sitional priors to get SPIRE fluxes (Roseboom et al. 2012) and then 
matched to the optically-based catalogues from Capak et al. (2007) 
and Ilbert et al. (2009). The same algorithm has been adopted 
for SPIRE-250 ^m sources with no PACS detection. Totally, the 
number of FIR sources in COSMOS with either PACS or SPIRE- 
250 pm detection is about 17000. 

Given that optical and near-IR observations in these fields are 
deeper than Far-IR ones obtained from Herschel, the fraction of 
//eri'c/ie/-selected galaxies without a counterpart in optical/near- 
infrared catalogues reaches only a few per cent in each field. 


2.2.1 X-ray counterparts 

We used optical/near-infrared counterpart positions to cross-match 
our Herschel-selected sample with available X-ray data from Chan¬ 
dra observations in COSMOS and in the GOODS fields. 

• In the GOODS-South, we use the 4-Ms Chandra-Deep Eield 
South (CDF-S) observations (Xue et al. 201 1). The X-ray catalogue 
provides count rates and observed fluxes for each source in different 
bands: soft (0.5-2 keV), hard (2-8 keV) and full (0.5-8 keV). 

• In the GOODS-North, X-ray data are available from the 2-Ms 
Chandra-Deep Field North (CDF-N) observations (Alexander et al. 
2003). 

• In the COSMOS field, observations from the Chandra- 
COSMOS (C-COSMOS, Elvis et al. 2009; Civano et al. 2012) sur¬ 
vey are publicly available but only cover about 0.9 deg^ (reaching 
about 160 ks in the central 0.45 deg^ and 80 ks outside) instead of 
~2 deg^ scanned by Herschel. Consequently, we limit the parent 
sample to match the common sky area, which implies a cut to ~45 
per cent of the original population, leading from roughly 17000 to 
7272 FIR-selected galaxies. 

After having collected as many spectroscopic redshifts as pos¬ 
sible (see Section 2.3), we cut our original sample at redshift z<2.5, 
since above this threshold poor statistics would affect the signifi¬ 
cance of our results. This leads to 829, 804 and 7011 sources in 
GS, GN and COSMOS, respectively, for a total number of 8644 
FIR-selected galaxies across all fields. More than 90 per cent of 
these FIR-selected galaxies are detected in at least two Herschel 
bands. The cross-match with X-ray detections has been made via 
a neighborhood algorithm, by assuming 1 arcsec matching radius 
between optical positions of the counterparts assigned to X-ray and 
Herschel sources, respectively. The number of X-ray detections is 
212/829 in GS, 134/804 in GN and 448/7011 in COSMOS. 

As highlighted in previous studies (e.g. Rosario et al. 2012), 
we confirm that the fraction of Herschel sources detected in X- 
rays is generally small and is a function of field, depending on the 
relative depth of X-ray/IR observations. The fraction of Herschel 
sources detected by Chandra is 26 per cent in the GOODS-S, 17 
per cent in the GOODS-N and only ~6 per cent in COSMOS. In 
all fields, we have taken the observed (i.e. obscured) X-ray fluxes 
in each band (soft, hard and full) from the publicly available cata¬ 
logues. 

We derived the average X-ray properties for the rest of the 
sample by performing a stacking analysis on X-ray maps (see Sec¬ 
tion 4.1 and Appendix A). 
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Figure 1. Redshift distribution of Herschel galaxies in COSMOS (top 
panel) and GOODS (bottom panel). Spectroscopic and photometric red- 
shifts are shown in red and blue, respectively, while the black line is the 
sum of the two. Note that the scale of the y-axis is logarithmic. 

2.3 Spectroscopic and photometric redshifts 

Extensive redshift compilations are publicly available in both 
GOODS and COSMOS. In Fig. 1 the redshift distributions for the 
GOODS and COSMOS samples are shown, distinguishing between 
spectroscopic (red) and photometric (blue) measurements. We de¬ 
fer the reader to Berta et al. (2011) for a careful description of red¬ 
shift catalogues, including uncertainties on photometric redshifts. 
Here we just provide a short list of references of redshift measure¬ 
ments. The redshift completeness for our Herschel sources with 
counterpart in optical/near-infrared catalogs reaches 100 per cent 
in all fields. 

• GOODS: in GOODS-S, we extended the original spectro¬ 
scopic sample presented by Grazian et al. (2006) and Santini et al. 
(2009) in the GOODS MUlti-wavelength Southern Infrared Cata¬ 
log (GOODS-MUSIC) with publicly available spectroscopic red¬ 
shifts, as described in Berta et al. (2011), reaching a global spec¬ 
troscopic fraction for our Herschel sample as high as 67 per cent. 
In GOODS-N, redshift measurements are taken from Berta et al. 
(2011), who collected spectroscopic redshifts from Barger et al. 
(2008) for about 64 per cent of the //erac/zeZ-selected sample. In 
both fields, photometric redshifts have been derived by Berta et al. 
(2011) by using the EAZY (Brammer et al. 2008) code, as de¬ 
scribed in Wuyts et al. (2011a). 

• COSMOS: we used photometric redshifts from Ilbert et al. 
(2010) and Wuyts et al. (2011a). For Chandra detected sources 
we have made a cross-check with photometric redshifts presented 
by Salvato et al. (2011) which are more suitable for AGN dom¬ 
inated sources. We retrieved spectroscopic measurements from 
the zCOSMOS survey by Lilly et al. (2007, 2009), either the 
public zCOSMOS-bright or the proprietary 2 :COSMOS-deep data 
base. We also browsed the most recent public spectroscopic sur¬ 
veys, replacing our photometric redshifts with spectroscopic ones 
in case of high reliability: Ahn et al. (2014), from Data Release 
10 of the Sloan Digital Sky Survey (SDSS); Coil et al. (2011) 


from the PRIsm MUlti-object Survey (PRIMUS) catalogue; Trump 
et al. (2009) from the COSMOS-Magellan spectroscopic catalogue. 
Globally, about 50 per cent of the COSMOS FIR-selected sample 
have a spectroscopic redshift. 

The overall fraction of spectroscopic redshifts is larger than 
60 per cent at z<1.5, while it decreases to 20 per cent at z~2. We 
stress that, if limiting our sample to spectroscopic redshifts only, 
all the results would remain consistent within la uncertainty with 
those already presented and discussed in Sections 5 and 6. 


3 THE SFR-M.PLANE 

Given the wealth of photometric datapoints available in all fields 
from the UV to the sub-mm, we performed broad-band SED de¬ 
composition to derive M* and SFR for the entire sample. Each ob¬ 
served SED has been fitted with the MAGPHYS^ code (da Cunha 
et al. 2008), as well as with a modified version of MAGPHYS 
adapted to include an AGN component (Berta et al. 2013) using 
AGN templates by Fritz et al. (2006) and Feltre et al. (2012). The 
best-fit obtained with the AGN is preferred if the resulting value 
significantly decreases (at >99 per cent confidence level, on the ba¬ 
sis of a Fisher test) compared to the fit without the AGN. In this 
case, SFRs and M* estimates are taken from the fit with AGN, 
otherwise from the original MAGPHYS code. However, even if the 
AGN component is required in the best-fit, we stress that its contri¬ 
bution to the galaxy IR luminosity is not dominant (i.e. lesssimlO 
per cent) for most of the sample. We defer the reader to Berta et al. 
(2013) for further details on SED decomposition, and to Delvec¬ 
chio et al. (2014) for statistical analysis. 

The SFR has been derived by converting the total IR (rest 8- 
1000 fj,m) luminosity taken from the best-fit galaxy SED (i.e. cor¬ 
rected for a possible AGN emission) using the conversion from 
Kennicutt (1998), scaled to a Chabrier (2003) IMF®. The M* is 
derived from the SED decomposition itself, which allows to get 
robust measurements also for type-1 AGN, as most of them show 
near-IR (~1 fim) emission dominated by the host galaxy light (e.g. 
Bongiorno et al. 2012). We have checked that our estimates of M* 
for optically identified type-1 AGN are consistent within the uncer¬ 
tainties (around 0.3 dex) with those presented by Bongiorno et al. 
(2012), with no systematics. 

The sample has been split in five different redshift bins: 
0.01<z<0.25, 0.25<z<0.50, 0.50<z<0.80, 0.80<z<L50 and 
L50<z<2.50. We place our galaxies on the SFR-M, plane in 
Fig. 2, marking with different colours GOODS (red) and COS¬ 
MOS (blue) sources. Given that PACS flux density limits in the 
GOODS fields are about 5 times lower than in COSMOS, Her¬ 
schel observations in the GOODS fields detect fainter IR galax¬ 
ies (i.e. lower SFR) compared to observations in COSMOS at the 
same redshift. Since our Herschel-based selection is sensitive to 
the most star-forming galaxies in each field, the wedge traced by 
the observed galaxy distribution in the SFR-M* is relatively flat 
(see Rodighiero et al. 2011) compared to the linear main-sequence 

MAGPHYS can be retrieved at http://www.iap.fr/magphys/ 
magphys /MAGPHYS . html 

® We computed the SFR for each galaxy by accounting for its obscured 
star formation only. As pointed out by Magnelli et al. (2013), the fraction 
of unobscured SFR density ranges between 12 per cent and 25 per cent at 
z<2, but drops to a few per cent for FIR-selected samples of galaxies (e.g. 
Wuyts et al. 201 la,b). 
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Figure 2. Individual fferic/ieZ-selected sources as a function of SFR, M» and redshift. Red and blue symbols refer to GOODS-S/N and COSMOS sources, 
respectively. The black solid line at any redshift represents the main-sequence defined by Elbaz et al. (2011). The black dashed line marks 4 X higher sSFR as 
a threshold between main and off-sequence galaxies. 


relation defined by Elbaz et al. (2011). The MS evolution with red¬ 
shift is parametrized as sSFRms = 26 x [Gyr“^] (Elbaz et al. 

2011), where tcosmic is the cosmic time (in Gyr) starting from the 
Big Bang. Our selection also includes off-sequence galaxies, with 
sSFR > TxsSFRms (~0.6 dex), corresponding to symbols above 
the black dashed lines. In each redshift slice, the sample has been 
split in both SFR and M*, taking 0.5x0.5 dex bins®. Instead of 
keeping a fixed binning configuration in SFR and M* at all red- 
shifts, we decided to center our bins on the main-sequence relation 
in each bin of M* and redshift^®. This arrangement is preferable, 
since it allows us to better highlight any potential systematics in 
terms of average AGN properties between main and off-sequence 
galaxies. 


4 X-RAY ANALYSIS 

In this Section we present the X-ray analysis of our Herschel se¬ 
lected sample. For sources detected in X-ray (~ 10 per cent) we 
get count rates, fluxes and observed luminosities in different bands 
from the above-mentioned catalogues (section 2.2.1). X-ray unde¬ 
tected sources represent most (~ 90 per cent) of the Herschel sam¬ 
ple studied in this work. To derive their average X-ray properties, 
we stack the X-ray maps (section 4.1), in the observed soft, hard 
and full bands. Subsequently, we characterize the average X-ray 
properties of Herschel galaxies in each bin of SFR, M* and redshift 
by considering X-ray detections and stacks together (section 4.2). 

® We note that the typical uncertainty on individual SFR and M* mea¬ 
surements is of the order of 0.2-0.3 dex, so significantly smaller than the 
bin-width. 

We checked also for alternative binning configurations (i.e. either dif¬ 
ferent shapes for individual bins, or different placement on the SFR-M* 
plane) but this did not produce any significant impact on our results. 


Then we subtracted the X-ray flux expected from star formation 
(section 4.3) in single X-ray bands and corrected the remaining X- 
ray emission for the nuclear obscuration (section 4.4) to derive the 
intrinsic (unobscured) mean AGN X-ray luminosity (rest-frame). 

Using widely-adopted conversion factors, we employ the final 
nuclear X-ray luminosity as a proxy to evaluate the mean BHAR in 
each bin of SFR, M* and redshift (see Section 5). 


4.1 Stacking the X-ray maps 

Here we briefly mention the main steps concerning the X-ray stack¬ 
ing and refer the reader to Appendix A for further details. We com¬ 
bined all Herschel sources undetected in X-rays from both GOODS 
and COSMOS fields and we grouped them as a function of SFR, 
M* and redshift. After masking all X-ray detected sources which 
could potentially affect the stacked signals, we piled up single 
cutouts of X-ray undetected sources in the same bin of SFR, M* 
and redshift, centered on their optical coordinates. For each object, 
we defined regions from which we extract source and background 
photons, in order to derive background-subtracted (i.e. net) pho¬ 
ton counts. Finally, we needed to correct for differential sensitiv¬ 
ities between various X-ray fields, so we normalized the resulting 
net photon counts of each object by the corresponding effective 
(i.e. corrected for instrumental effects) exposure time, which pro¬ 
vides exposure-corrected mean count rates in each observed X-ray 
band^^. To convert count rates into observed fluxes, we assumed a 
power-law spectrum with r= 1.4, based on the empirically observed 
spectrum of the X-ray background (e.g. Gilli et al. 2007). 


We note that our results are in good agreement, within the uncer¬ 
tainties, with those obtained using CSTACK (http: //cstack . ucsd . 
edu/cstack/, developed by T. Miyaji) on the same X-ray maps. 
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4.2 Mean X-ray luminosity of Herschel sources 

To get the overall X-ray properties of our Herschel-selected galax¬ 
ies, we collected in each bin of SFR, M* and redshift both X-ray 
detected and undetected sources. Assuming that the bin includes 
N sources, m detected and n undetected in X-rays, we computed 
a number weighted average of their X-ray fluxes, according to the 
following formula: 


L^"‘^[erg s"^] = 5.13 x 10^® SFR [M© yr"^] (3) 

In each bin we subtracted the X-ray luminosity due to star forma¬ 
tion from the average X-ray luminosity (section 4.2), both in the 
soft and the hard band. However, we note that this subtraction does 
not impact significantly the previous values, as the typical non- 
AGN contribution arising from the joint (X-ray detected -l- stacked) 
sample is less than 10 per cent. 


(S'bin) 


n X Satack + ‘S'i 

N 


( 1 ) 


where 5'atack is the mean X-ray flux obtained by stacking n un¬ 
detected sources and Si is the individual flux measured for the i- 
th detected source within the same bin. We stress that the above- 
mentioned expression is a number weighted average, therefore ap¬ 
propriate to investigate mean properties (e.g. mean X-ray luminos¬ 
ity) of the underlying galaxy population, but we caution that it does 
not necessarily represent the most probable value of the observed 
distribution of a given parameter. Considering that ~90 per cent of 
Herschel galaxies are not X-ray detected, as well as the fact that X- 
ray detections are about 50 times brighter than stacks, in most bins 
X-ray detected and stacked sources provide comparable contribu¬ 
tions to the mean X-ray flux. Rest-frame average X-ray luminosi¬ 
ties in each band, for both X-ray detected and stacked sources, are 
derived from mean X-ray fluxes by assuming a power-law spectmm 
with intrinsic slope r=1.9 (Tozzi et al. 2006; Mainieri et al. 2007) 
and no obscuration. We performed this calculation as a function of 
SFR, M* and redshift. 


4.3 Subtraction of X-ray emission from star formation 

Using X-rays to investigate the level of AGN activity in the SFR- 
M* plane requires the subtraction of X-rays from other processes 
in the host galaxy, especially those related to star formation. This 
subtraction will allow us to apply the correction for nuclear obscu¬ 
ration (Section 4.4) only to the AGN-related X-ray emission rather 
than the total (i.e. AGN -l- galaxy) one. It is known that X-ray ob¬ 
servations provide a relatively clean selection of accreting SMBHs, 
but a fraction of the X-ray emission (e.g. in Chandra Deep Fields, 
CDFs) may also come from X-ray binaries and the hot interstellar 
medium (e.g. Mineo et al. 2012a,b). A single threshold (i.e. Lx = 
3 X10'*® erg s“^) in X-ray luminosity is not a proper way of classi¬ 
fying these galaxy populations, since it is known that X-ray emis¬ 
sion from star-forming galaxies shows a positive correlation with 
SFR. Previous works (e.g. Ranalli et al. 2003) calibrated this rela¬ 
tion in the local Universe, while recent studies (Vattakunnel et al. 
2012; Mineo et al. 2014; Symeonidis et al. 2014) uncovered this 
relation up to z~1.5. 

Though all these relations are calibrated through independent 
analyses and selection techniques, they provide reasonably consis¬ 
tent (within a factor of 2) estimates of the X-ray emission expected 
from star formation. However, we prefer the relation by Symeonidis 
et al. (2014), since they also exploited Herschel data and performed 
stacking on X-ray maps to better characterize the average Lx-SFR 
correlation in inactive (i.e. non-AGN) SFR-selected galaxies. They 
found a quasi linear relation holding at 1<SFR<1000 Mq yr“* 
and not evolving significantly with redshift up to z~1.5. We scaled 
their relation to a Chabrier (2003) IMF and converted to our rest- 
frame soft (0.5-2 ke'V) and hard (2-8 keV) X-ray bands: 

L^°‘'‘[erg s"^] = 2.04 x 10®® SFR [Mq yr"*] (2) 


4.4 Correction for nuclear obscuration 

To derive the intrinsic AGN luminosity in each bin of SFR, M* and 
redshift, we considered the AGN-related X-ray emission derived 
in Section 4.3 and followed the approach developed by Xue et al. 
(2010), who used the hardness ratio (HR) as a proxy to estimate the 
nuclear obscuration. The hardness ratio is defined as: 

HR = (CRhard “ CRaoft)/ (CRhard + CRgoft) (4) 

where CRhard and CRsoft represent the (exposure-corrected) 
count rates in hard and soft-bands, respectively. 

Given that we are not able to constrain the HR for each in¬ 
dividual galaxy of the sample, we correct the mean AGN X-ray 
emission of all Herschel galaxies (both detected and undetected in 
X-rays, see Eq. 1) for an average level of obscuration. We make the 
simple assumption that the mean HR calculated after the subtrac¬ 
tion in soft and hard X-ray bands (Eqs. 2, 3) is representative of the 
galaxy population in the same bin of SFR, M* and redshift. 

We parametrized the effect of nuclear obscuration by assum¬ 
ing a single power-law X-ray spectrum, with intrinsic photon index 
F = 1.9 (model waxzwaxpo in XSPEC, Arnaud 1996) and ac¬ 
counting for absorption, both Galactic and intrinsic to the AGN. 
However, since the hardness ratio deals with photon counts instead 
of fluxes, one needs to convolve the intrinsic model with the instru¬ 
ment response curve.*® After performing this convolution, we built 
a set of simulated spectra that predict the observed hardness ratio as 
a function of hydrogen column density (Nh) and redshift. We cal¬ 
culated the observed hardness ratio from AGN count rates obtained 
in Section 4.3 and selected the spectral model (i.e. intrinsic column 
density) best reproducing the observed HR at the mean redshift of 
the underlying galaxy population in the same bin of SFR and M*. 
From the observed-frame, absorption-corrected fluxes in the 0.5-8 
keV band (^[o.s-gi.int). we calculated the rest-frame, intrinsic full- 
band X-ray luminosity as follows: 

L[0.5-S],int = TtT D® S[o.5-8],int (1 + z)'" ® (5) 

where Dl is the luminosity distance corresponding to the mean red¬ 
shift of the Herschel population in a given bin, while the intrinsic 
photon index F is set to 1.9. 

As a sanity check, we also compared our obscuration- 
corrected full-band (0.5-8 keV) luminosities with those presented 
by Xue et al. (2011) in the Chandra-Deep field South (CDF-S) 4- 
Ms catalogue and found an excellent agreement, as expected, given 
that we followed similar approaches. We found that obscuration 
level does not significantly affect the average AGN X-ray luminos¬ 
ity, as the typical correction factor is about 1.3. The obscuration- 
corrected X-ray luminosities estimated through hardness ratio are 


The instrument response curve has been con'ected for several instru¬ 
mental effects: vignetting, exposure time variations, energy-dependent effi¬ 
ciency, variation of effective area with time of observations. 
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generally consistent within a factor of ~30 per cent with more pre¬ 
cise measurements from spectral-fitting analysis (Xue et al. 201 1). 
Nevertheless, as argued by Xue et al. (2011), the level of obscu¬ 
ration might be affected by strong uncertainties in case of highly 
obscured AGN. In addition, we compared our predictions based on 
HR with spectral measurements of ~400 Chandra-COSMOS AGN 
presented by Lanzuisi et al. (2013). We ended up with reasonably 
small scatter (about 0.2 dex) in [0.5-8] keV X-ray luminosities, 
even for sources classified as highly oscured (Nh > 10^^ cm~^) 
AGN from X-ray spectral-fitting. Given these sanity checks and 
the relatively small average corrections for obscuration obtained 
for our sample, this might suggest that the average intrinsic X-ray 
luminosity of Herschel-selected galaxies does not arise primarily 
from highly obscured AGN. However, a thorough X-ray spectral 
analysis of these sources would be beyond the scope of this paper. 


4.4.1 Uncertainty on the [0.5-8] keV intrinsic Lx 

We evaluated the uncertainty on the [0.5-8] keV intrinsic X-ray lu¬ 
minosity L[o, 5 _ 8 ],int by performing a bootstrapping analysis. This 
technique provides reliable error bars, especially in case a small 
fraction of the objects dominate the signal. 

Suppose there are N objects populating a given bin of SFR, 
M* and redshift, out of which m are detected in X-rays, while n 
are not detected. We selected at random N objects from the same 
bin, allowing duplication of the same source. This random extrac¬ 
tion likely leads to different numbers of detections (m') and non¬ 
detections (n'). We combined the photon counts together from the 
new n' sources to derive stacked count rates and fluxes. Then we 
applied Eq. 1 (see Section 4.2) to get the average X-ray flux rep¬ 
resentative of this random realization. With 1000 iterations of this 
calculation, we obtained the distribution of average count rates and 
fluxes of the galaxy sample. The 16*^ and 84*^ percentiles of the fi¬ 
nal distribution set the la lower and upper bounds on the measured 
X-ray flux. 

To evaluate the uncertainty on the obscuration-corrected Lx, 
we iterated 1000 times the same analysis described in Sections 4.3 
and 4.4, once for each random realization. This approach returns 
the final distribution of the [0.5-8] keV intrinsic X-ray luminosi¬ 
ties, with ±1 ct error bars estimated through bootstrapping. Since 
this technique would provide quite large error bars (i.e. almost un¬ 
constrained fluxes and luminosities) in case of poor statistics, we 
have required a minimum number of sources in each bin (> 15, 
regardless of the number of detections and non-detections). 


5 RESULTS 

We used the intrinsic X-ray luminosity as a proxy for mapping the 
average BHAR in the SFR-M, plane and studying its correlation 
with integrated galaxy properties. 


5.1 Average BHAR 

Obscuration-corrected X-ray luminosities in each bin have been 
turned into average bolometric AGN luminosities (Lboi) by as¬ 
suming a set of luminosity-dependent X-ray bolometric correc¬ 


tions from Marconi et al. (2004)^^. We have assumed a con¬ 
stant conversion factor to convert the average AGN bolometric 
luminosity (Lboi) (in erg s“^) to average SMBH accretion rate 
(Afbh(Mt, SFR, z)) (in MQyr~^), according to the formula (see 
Alexander & Hickox 2012): 


(Mbh(M., SFR, z)) =0.15 ^ 


(Lboi(M.,SFR,z)) 

10 « 


( 6 ) 


where the matter-to-radiation conversion efficiency e is assumed to 
be 10 per cent (e.g. Marconi et al. 2004). Table B1 lists redshift, 
SFR, M* and accretion rates for all bins involved, as well as the 
number of sources entering each bin. In Fig. 3 we show the SFR- 
M* plane at different redshifts with colour-coded average BHAR. 
A few upward and downward triangles set Icr lower and upper lim¬ 
its, respectively. They replace the formal values in case the correc¬ 
tion for obscuration is not applicable. Indeed, the subtraction of the 
star-formation X-ray emission from the total (i.e. AGN and star- 
formation) X-ray luminosity obtained in Section 4.2, in some bins 
left solely the hard (soft) X-ray emission, which returned a lower 
(upper) limit in the [0.5-8] keV intrinsic X-ray luminosity. 

In addition, an upper limit is also imposed in case the observed 
X-ray spectrum is softer than any spectral model with intrinsic 
r=L9^^. Indeed, this slope represents just a mean value among the 
overall distribution of AGN X-ray spectra, whose typical scatter is 
around 0.2 (e.g. Piconcelli et al. 2005). We note that our SFR-M*-z 
grid follows the evolution with redshift of the MS. Even at our large 
source statistics, uncertainties in BHAR remain noticeable (about 
0.4 dex on average, see Table B 1), also due to the fact that in some 
bins a few, intrinsically luminous X-ray detections might dominate 
the global signal, as stated in Equation 1 . 

In Fig. 3, the high SFR Herschel galaxies at higher redshift 
show larger average BHAR than the more local and lower SFR 
galaxies. Also, within the panels for each redshift slice, trends are 
indicated with SFR and/or M,. We proceed to study which of these 
galaxy properties correlates best with BHAR, and discuss results in 
the context of 0<z<2.5 galaxy evolution. 


5.2 Correlation of BHAR with galaxy properties 

We explore here the observed trends between (Mbh) and various 
galaxy properties: SFR, M* and offset from the MS^®. Because of 
the evolution with cosmic time of SFRs and BHARs, we fit the data 
separately for each redshift bin, by assuming a linear trend in the 
log-log space through: 

log((Mbh)) = a X log(a;) -I- /3 (7) 

where a and 13 represent slope and intercept, respectively, while 
X is the corresponding independent variable. We have used the 
IDL routine lMSL_MULTlREGRESS.PRO, which performs a linear 
regression fit considering error bars in both variables, and returns 
best-fit intercept and slope with related Icr uncertainties. The linear 
best-fits obtained in each redshift bin are summarized in Table 1 . 

We remark that if taking a fixed bolometric correction value of 22.4 (as 
done by Mullaney et al. 2012a and Chen et al. 2013), the AGN bolometric 
luminosities would be larger by a factor of about 2, for Lx < lO'^'^ erg s~ ^. 

The softest X-ray spectrum that is reproducible with intrinsic slope 
r=1.9 corresponds to HR = -0.47. 

As the adopted MS relation shows a linear trend at all redshifts, the 
offset from the MS becomes simply the ratio between the sSFR and that 
corresponding to the main-sequence sSFRms. for a given bin of M* and 
redshift. 
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Figure 3. BHAR distribution (colour-coded) in the SFR-M* plane at 0.01<z<2.5. The bins in SFR and M* are arranged to sample the main-sequence 
relation (black solid line), which evolves with redshift. The black dashed line divides main-sequence from off-sequence galaxies. Coloured bins include at 
least 15 sources, either detected or undetected in X-rays. Upward and downward triangles set \a lower and upper limits, respectively, on the average BHAR. 
In particular, as the expected non-AGN contribution was subtracted from the average X-ray luminosities obtained in Section 4.2, in some bins this subtraction 
left solely the hard (soft) X-ray emission, which returned a lower (upper) limit in the final intrinsic X-ray luminosity. 


Table 1. Best-fit parameters returned by fitting a linear relation in the log-log space, as a function of redshift, between {Mhh) ^^d galaxy properties {x- 
parameter). a and /3 are the slope and intercept of the linear best-fit (see Eq. 7). Values in brackets set the Icr uncertainty on the related parameters. The 
Spearman’s rank correlation coefficient p is the strength of the correlation, while P(correlation) represents the significance of its deviation from zero, consid¬ 
ering the number of points N in each panel. In columns 5 and 6 their estimates have been derived from a linear regression fit, while in columns 7 and 8 we 
used a partial-correlation fitting (see text for details). 


(1) j 2 -bin 

(2) x-parameter 

(3) a 

(4)/3 

(5)P 

(6) P(correlation) 

(7) partial p 

(8) partial P(correlation) 

{9)N 

0.01 < z < 0.25 

SFR [Mo yr-l] 

0.54 (±0.27) 

^.31 (±0.12) 

0.68 

95.8% 

0.39 

83.0% 

9 

0.25 < z < 0.50 

SFR[M0 yr-l] 

1.06 (±0.24) 

-4.16 (±0.38) 

0.75 

99.1% 

0.49 

92.6% 

11 

0.50 < z < 0.80 

SFR[M0 yr-i] 

1.05 (±0.52) 

^.02 (±0.64) 

0.60 

96.1% 

0.40 

89.9% 

12 

0.80<z< 1.50 

SFR [Mo yr-i] 

1.44 (±0.30) 

^.29 (±0.46) 

0.80 

99.9% 

0.71 

99.6% 

14 

1.50<z<2.50 

SFR [Mo yr-l] 

1.13 (±0.38) 

-3.42 (±0.78) 

0.68 

99.0% 

0.68 

99.2% 

13 

0.01 < z < 0.25 

M. [Mo] 

0.44 (±0.21) 

-8.63 (±2.17) 

0.73 

97.5% 

0.52 

91.7% 

9 

0.25 < z < 0.50 

M. [Mol 

0.76 (±0.33) 

-11.41 (±3.51) 

0.75 

99.1% 

0.52 

93.8% 

11 

0.50 < z < 0.80 

M. [Mol 

1.11 (±0.34) 

-14.59 (±3.66) 

0.58 

95.2% 

0.41 

89.9% 

12 

0.80<z< 1.50 

M. [Mol 

1.07 (±0.18) 

-13.58 (±1.96) 

0.65 

98.8% 

0.43 

93.0% 

14 

1.50<z<2.50 

M. [Mol 

0.25 (±0.27) 

-3.89 (±2.94) 

0.23 

56.3% 

-0.17 

30.0% 

13 

0.01 < z < 0.25 

sSFR/sSFRnis 

-0.09 (±0.27) 

^.17 (±0.13) 

-0.02 

33.2% 

/ 

/ 

9 

0.25 < z < 0.50 

sSFR/sSFRnis 

-0.18 (±0.42) 

-3.33 (±0.21) 

-0.17 

43.6% 

/ 

/ 

11 

0.50 < z < 0.80 

sSFR/sSFRms 

-0.53 (±0.50) 

-2.74 (±0.17) 

0.03 

36.7% 

/ 

/ 

12 

0.80<z< 1.50 

sSFR/sSFRms 

-0.53 (±0.40) 

-2.13 (±0.15) 

0.19 

54.6% 

/ 

/ 

14 

1.50<z<2.50 

sSFR/sSFRnis 

0.10 (±0.31) 

-1.10 (±0.15) 

0.40 

83.8% 

/ 

/ 

13 


5.2.1 Spearman’s rank 

To evaluate the significance of the observed trends, we used 
the Spearman’s rank correlation coefficient p, that indicates the 
strength of any observed correlation. However, the Spearman’s test 


considers the observed data points as “exact” and does not take 
into account their possible error bars. To estimate the most proba¬ 
ble correlation coefficient given the error bars, we used the proce¬ 
dure detailed in Santos-Sanz et al. (2012) (see their Appendix B.2 
for details). Briefly, they generated 1000 samples of data points, 
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Figure 4. Average BHAR vs SFR, at different redshifts. The colour-coded 
bar is the M*, while individual black dashed-dotted lines show the best-fit 
linear relation inferred in each z-bin, whose coefficients are listed in Table 
1. Bins representing off-sequence galaxies are enclosed in black squares. 
We stress that the lowest redshift bin might suffer from incompleteness at 
high SFRs. 


Figure 5. Average BFIAR vs M*, at different redshifts. The colour-coded 
bar is the SFR, while individual black dashed-dotted lines set the best-fit 
linear relation inferred in each z-bin, whose coefficients are listed in Table 
1. Bins representing off-sequence galaxies are enclosed in black squares. 
We stress that the lowest redshift bin might suffer from incompleteness at 
high M*. 


building each synthetic dataset from its associated Gaussian distri¬ 
bution function, where error bars correspond to one standard devi¬ 
ation. This Monte Carlo bootstrap analysis provides a distribution 
of Spearman’s rank correlation coefficients, where its median value 
gives the weighted Spearman’s coefficient that represents the most 
likely correlation to the observed data points. We used the IDL rou¬ 
tine R_CORRELATE.PRO to derive the average value of Spearman’s 
p: in addition, this function provides the significance of its devia¬ 
tion from zero P (correlation) of the observed trend, considering 
the number of points populating the correlation. These parameters 
are used to compare different relationships and infer the strength of 
the observed correlations. We applied this analysis when studying 
the relationship between average BHAR and SFR, M* and MS- 
offset in each redshift bin (see Table 1). 

5.2.2 Correlation of BHAR with SFR and M, 

Fig. 4 and 5 show the relations (Mbh) vs SFR and (Mbh) vs M*, 
respectively, in different redshift bins. Circles correspond to the 
same bins already shown in Fig. 3, but projected on the (Mbh)- 
SFR or the {Mbh)-M, plane. Error bars of each bin mark the Tier 
uncertainties in (Mbh), obtained by applying a Monte Carlo boot¬ 
strapping (see Section 4.4.1) to the [0.5-8] keV intrinsic X-ray lu¬ 
minosity and rescaling to (Mbh) as described in Section 5.1. The 


black dashed dotted line in each z-bin represents the best-fit linear 
relation parametrized in Eq. 7, whose parameters are listed in Ta¬ 
ble 1. The best-fit slope (a) suggests that the correlation between 
BHAR and SER is consistent with a linear relation, at z>0. 25. 
Eig. 4 and 5 show that the best-fit relations at 0.01<z<0.25 are 
flatter than linear. This is likely due to the fact that the comov¬ 
ing volume covered by our fields is not large enough to detect Ul¬ 
tra Luminous InfraRed Galaxies (ULIRGs, Lm >10 ^^Lq) hosting 
powerful Quasars. At z>0.25 both relationships show clear trends 
with nearly linear slopes and relatively high significance (>2a 
level). Eocusing on individual redshift slices, the comparison be¬ 
tween Spearman’s rank coefficients in Table 1 (columns 5 and 6) 
shows that the strength of the observed correlations of BHAR with 
either SER or M* are comparable at z<0.8, while at higher redshifts 
the correlation with SER becomes progressively more significant. 
Especially at z~2, these coefficients significantly favour the corre¬ 
lation with SER, and show a nearly flat trend between (Mbh) and 
M*. 

Another way of illustrating this point is by looking at the 
z>0.5 panels in Eig. 4 and 5; Eig. 5 shows that off-sequence galax¬ 
ies (highlighted with black squares) have higher average BHAR, 
compared to mass-matched main-sequence galaxies. In contrast, at 
a given SER in Eig. 4, the BHAR for off-sequence galaxies does not 
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Table 2. Best-fit parameters obtained in each z-bin through multiple- 
correlation fitting procedure. 7 is the slope of SFR, while <5 refers to M*. 
Numbers in brackets represent their Icr uncertainties. The last column pro¬ 
vides with the reduced values returned from the multiple linear fitting. 


2 -bin 

7 (SFR) 

(S(M*) 

xl 

0.01 < z < 0.25 

0.38 (±0.20) 

0.31 (±0.16) 

1.90 

0.25 < z < 0.50 

0.78 (±0.32) 

0.52 (±0.24) 

1.78 

0.50 < z < 0.80 

0.48 (±0.36) 

0.92 (±0.29) 

1.93 

0.80 < z < 1.50 

0.85 (±0.22) 

0.73 (±0.16) 

0.89 

1.50 < z < 2.50 

1.16 (±0.33) 

-0.04 (±0.18) 

1.86 


stand out of those for main-sequence galaxies. Overall, this may in¬ 
dicate that at z>0.5 the BHAR is primarily dependent on SFR. 

Given the well known main-sequence relation between SFR 
and M* at all redshifts considered here, it is necessary to carry this 
analysis further and simultaneously investigate the dependence of 
{Mhh} on both SFR and M*. We addressed this issue by perform¬ 
ing a partial-correlation fitting analysis. We used the IDL routine 
R_CORRELATE.PRO to evaluate the Spearman’s p related to each 
couple of parameters and then we combined them according to the 
following expression: 


Pabc — 


Pab Pacpbc 

7Ti - PicTa^’-TO 


( 8 ) 


which returns the partial correlation between a and b adjusted for 
c. Assuming that {a,b,c) = (SFR, (Mhh), M*), then pabc would 
represent the partial correlation between SFR and {Mhh}, remov¬ 
ing the concomitant dependence on M*. As done in Section 5.2.1, 
we performed again a Monte Carlo bootstrap analysis to provide 
a distribution of Spearman’s coefficients from partial-correlation 
analysis, taking the median value of that distribution as the most 
probable to reproduce our data'^®. The latter p values obtained from 
partial-correlation fitting (herafter “partial p”), together with their 
corresponding significance levels, are listed in Table 1 (columns 7 
and 8 ). 

The comparison between partial p values and those obtained 
from a rank correlation of individual parameters shows similar 
trends. This suggests that the main-sequence relation between SFR 
and M* does not significantly affect the correlations with the aver¬ 
age BHAR. The observed trends show that a correlation with M* 
is slightly preferable than with SFR at z<0.5, while at higher red- 
shifts the correlation with SFR becomes progressively more sig¬ 
nificant. However, we evaluated the difference between SFR and 
M*, in terms of resulting partial-correlation coefficients, to be less 
than la at z<0.8 and of the order of 2.2a at 1.5<z<2.5. This sug¬ 
gests that for our current sampling of the SFR-M* plane, the over¬ 
all evolution of the average BHAR is best represented through a 
joint dependence on both SFR and M, at all redshifts. A difference 
is detected only in the highest z-bin, where the trend with SFR is 
fairly preferable with respect to that with M*. 

We performed a multiple-correlation fitting in the log space 
to provide a simple analytic expression of the simultaneous depen¬ 
dence of {Mhh) on both SFR and M* as a function of redshift. 
According to the following parametrization: 


log((Mbh)) oc 7 log(SFR) -f 51og(M.) (9) 


If taking the mean value instead of the median would change the re¬ 
sulting p by around 2-3 per cent. Moreover, we note that, if the observed 
trend is poorly significant, the most probable p value might prefer a mildly 
negative rather than a mildly positive correlation (see Table 1), or viceversa. 
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Figure 6. Average BHAR vs offset from the MS, at different redshifts. The 
colour-coded bar is the SFR, while individual black dashed-dotted lines set 
the best-fit linear relation inferred in each z-bin, whose coefficients are 
listed in Table 1. Bins representing off-sequence galaxies are enclosed in 
black squares. 


we defined 7 and 5 as the slopes of SFR and M*, respectively, and 
listed their values in Table 2 for each z-bin. At z<0.8, the slopes re¬ 
lated to SFR and M* are always consistent with each other (within 
la uncertainty), meaning that both parameters are required to trace 
the evolution of the average BHAR. As expected from our previous 
analyses, at z ~2 the correlation with {Mhh) is primarily driven by 
the SFR. 

Given that we are studying Herschel-selected galaxies, our 
sample is biased towards the most star-forming galaxies at any 
given redshift. Since this selection effect could potentially bias our 
results, we tried to evaluate the implications that incompleteness 
effects might have on this study. We repeated the same analysis by 
limiting our parent sample to Herschel galaxies with FIR flux (in 
either BAGS or SPIRE bands) larger than the flux corresponding 
to 80 per cent completeness level. This cut strongly reduces our 
statistics (about a factor of three), as well as the range of galaxy 
properties we are allowed to probe, but it allows us to test the valid¬ 
ity of our findings in a reasonably complete sample of star-forming 
galaxies. We repeated the same analysis and evaluated the strength 
of the observed correlations between average BHAR, SFR and M*, 
on the basis of this quite complete sample of star-forming galaxies. 
We ended up with consistent trends (within Icr uncertainty) in all 
redshift bins, with overall correlations between {Mhh) and galaxy 
properties in agreement with our previous analysis. In addition, we 
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checked that the observed trends obtained separately in GOODS 
fields and in COSMOS provide consistent results. 

We have found in this Section that average BHAR correlates 
with SFR and with M* for star-forming galaxies in a wide range 
of SFRs and redshifts. SFR may be the primary driver of these 
trends, at least for z>0.8. A further investigation will be feasible 
through forthcoming Chandra observations of the entire COSMOS 
field from the COSMOS Legacy Survey (PI: F. Civano), by which 
we plan to increase statistics and significance of the correlations in 
future work. 

5.2.3 Correlation of BHAR with offset from the MS 

Given long-standing evidence about AGN in local massive above- 
MS objects (ULIRGs), it is natural to probe for a link of BHAR 
and offset from the MS. This is shown in Fig. 6 and the best-fit 
parameters obtained through a linear regression fitting are listed 
in Table 1. It is evident that the average BHAR does not depend 
significantly on the offset from the MS. At all redshifts, the best- 
fit slope is consistent with a fiat or slightly declining relation. Our 
analysis based on Spearman’s p coefficients suggests a very weak, 
or absent, correlation at any redshift here considered^^. 

If we focus our analysis only to the few bins corresponding 
to above-MS objects (black squares in Fig. 5), they tend to show, 
at a given M*, larger BHAR compared to galaxies placed on the 
MS relation at the same redshift. However, this increase of BHAR 
from main-sequence to off-sequence galaxies is not observed when 
considering all bins, probably due to the wider range in specific 
SFR covered by MS galaxies in our sample. In addition, we have 
looked into possible trends of the ratio BHAR/SFR with offset from 
the MS and found no significant one at any redshift. 

In Section 6 we discuss this lack of trends in the framework 
of classical merger scenarios. All these findings are consistent with 
a link of BHAR to SFR, irrespective of whether a given SFR is 
reached in a bursting above-MS galaxy of moderate M,, or a high 
M* main sequence galaxy. 

5.3 BHAR as a function of SFR and M* over 0.01<z<2.5 

After analysing redshift bins separately, we now merge results over 
the full redshift range covered by our study. The top panel of Fig. 7 
shows the average BHAR (combining X-ray detected and stacked 
sources) in bins of SFR and M, from all our redshift slices. We have 
restricted ourselves here to the subsample with >80 per cent com¬ 
pleteness in SFR. The central and bottom panels show the projec¬ 
tions of the previous plot on the SFR-(Mbh) plane and M*-(Mbh) 
plane, respectively, with a redshift colour-coding. 

In these panels, the data from different redshifts better connect 
into a tighter and more consistent sequence for BHAR=f(SFR) than 
for BHAR=f(M*). This is in support of a primary dependence of 
the BHAR on SFR, though noting the less pronounced difference 
seen above within individual redshift slices. 

5.4 Comparison with previous studies 

We compare the relationships between AGN accretion and galaxy 
properties obtained from our analysis with those found by Mul- 

We have verified that our results remain if using the MS definition of 
Whitaker et al. (2012) or Schreiber et al. (2014), instead of the one proposed 
by Elbaz et al. (201 1). 
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Figure 7. Top panel: distribution of bins including Herschel galaxies above 
the >80 per cent completeness limit, as a function of BHAR, SFR and M* 
at all redshifts. The x, y axes represent the SFR-M* plane, while the 2:-axis 
sets the colour-coded BHAR. Note that the three-dimensional space is log¬ 
arithmic. Middle panel: previous plot projected on the SFR-(Mbh) plane, 
with colour-coded redshift. Bottom panel: projection on the M*—(Mth) 
plane, with colour-coded redshift. 


laney et al. (2012a, M12 hereafter) and Chen et al. (2013, C13 
hereafter). M12 analysed a sample of star-forming galaxies in the 
GOODS-S using BzK and 24 pm selection. The authors claimed a 
linear relation at z~2 between (Mbh) and M*, while Fig. 5 shows 
a weak trend. We checked whether the lack of trend in our data 
may be due to our sample selection. We match the selection crite- 
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ria adopted by M12, by including main-sequence galaxies only (i.e. 
sSFR < 4x sSFRms), and splitting the sample at 1.5<z<2.5 only in 
bins of M*. In addition, we used a constant bolometric correction of 
22.4 (from Vasudevan & Fabian 2007) to derive the average BHAR. 
We found that (Mbh) oc while M12 derived (Mbh) oc 

M105±o. 36 j^Q^yever, this nearly linear relation vanishes as we 
split the sample as a function of both SFR and M*: keeping the 
same Ml 2-like sample and changing only the grid configuration, 
the resulting best-fit relation with M* is (Mbh) oc 

To study the simultaneous dependence of (Mbh) on both SFR 
and M*, we performed a multiple-correlation fitting in the log 
space, as a function of (Mbh), SFR and M* at the same time. We 
used a Ml 2-like sample, but splitting it as a function of both SFR 
and M*. We found that (Mbh) oc These 

slopes may suggest a “pure” SFR-driven relation. 

In summary, differences between our findings and M12 are 
due to the different binning and analysis: for an MS-only mass- 
binned sample where M* and SFR are degenerate because of the 
MS relation, we reproduce their close to linear relation of M, and 
BHAR, while covering the SFR-M* plane in both dimensions al¬ 
lows to separate the two variables and emphasizes the role of SFR. 

Rodighiero et al. (2014, in preparation) analyse a mass- 
complete sample of galaxies at z~2 in the COSMOS field. We refer 
the reader to their work for a detailed comparison with results pre¬ 
sented by M12. 

We made also a comparison with results from C13, who 
measured a nearly linear trend between (Mbh) and SFR at 
0.25<z<0.80, splitting their Herschel sample in bins of SFR only. 
C13 found log BHAR = (-3.72±0.52) -l- (1.05±0.33) log SFR (both 
in Mq yr~^). By adjusting our grid configuration in bins of SFR 
to match their analysis, our fitting routine gives consistent coef¬ 
ficients, with BHAR = (-3.65±0.12) -l- (1.18±0.11) log SFR. We 
note that a similar trend is found even in case we split the sample 
in bins of SFR and M*. 


6 DISCUSSION 


In Section 5.2, we described the relationship between BHAR and 
SFR, M* and offset from the main-sequence, over a wide range 
of galaxy properties and cosmic epochs, from z~0 up to z=2.5. In 
this Section we discuss and interpret our findings in the context of 
current AGN/galaxy evolutionary scenarios. 

Mapping the average AGN accretion along both SFR and M* 
allows us to separate the degeneracy between the two parameters, 
shown by the presence of the MS relation since z~2. In Section 
5.2.2 we found that the average BHAR depends on both SFR and 
M* with similar significance levels. A more clear correlation with 
SFR is observed at z>0.8 and becomes more significant at z~2. 
This finding may support the SFR as the original driver of the cor¬ 
relation with average AGN accretion. 

We relate the trends of BHAR with SFR, M* and redshift to 
gas content and the redshift evolution of the gas fraction: 


/gas — 




Mg-as + M* 


( 10 ) 


which represents the ratio between total (i.e. molecular -l- neu¬ 
tral) gas mass and total baryonic mass of a galaxy. Several studies 


This trend is similar to what would be obtained in Fig. 6 (highest red¬ 
shift bin) by removing the two bins corresponding to off-sequence galaxies 
(black squares). 


have investigated the evolution of /gas in main-sequence galaxies 
from the local Universe to z~l (Leroy et al. 2008; Geach et al. 
2011; Saintonge et al. 2011a, 2011b, 2012), at z~1.5-2 (Daddi 
et al. 2010a; Tacconi et al. 2010, 2013), finding a strong (/gas oc 
(1-l-z)^) evolution up to z~2 and a plateau at higher redshift (z~3, 
Magdis et al. 2013). Direct CO observations presented by Daddi 
et al. (2010a) and Tacconi et al. (2010, 2013) provided first em¬ 
pirical evidence for the existence of very large gas fractions in 
z~l-2 main-sequence galaxies, with mass in gas even larger than 
the mass in stars (i.e. /gas >0.5). These gas rich systems have 
also higher SFRs (Daddi et al. 2010b; citealtGenzel-l-10), as ex¬ 
pected from the Schmidt-Kennicutt relation (Schmidt 1959; Ken- 
nicutt 1998). This gas-dominated phase in z~l-2 main-sequence 
galaxies is also reflected in a difference in morphology. Indeed, 
while local main-sequence galaxies are preferentially regular disks, 
their z~2 analogs show a larger fraction of irregular morphologies 
and/or clumpy disks (e.g. Elmegreen et al. 2007; Forster Schreiber 
et al. 2009; Kocevski et al. 2012), which make them potentially 
more efficient in funnelling the cold gas inward onto the SMBH 
(Bournaud et al. 201 1). 

In the light of the above-mentioned literature, it is justified to 
expect that the relationship between average BHAR and SFR be¬ 
comes stronger with increasing redshift, where the fraction of ac¬ 
tively star-forming gas becomes larger. Moreover, given that z~l-2 
galaxies are truly gas rich, especially those at lower M* (e.g. San- 
tini et al. 2014), it is also reasonable to assume that in these sys¬ 
tems the fraction of primordial gas (i.e. not yet converted into stars, 
therefore not causally linked to the galaxy stellar mass) is larger 
compared to local galaxies and to higher stellar mass galaxies at 
the same epoch. This may explain the weaker correlation of BHAR 
with M* that we found at higher redshift. 

In line with these arguments, Vito et al. (2014) found that 
AGN hosts are significantly more gas rich than inactive galaxies, 
at a given M* and redshift, suggesting that the probability that a 
SMBH is active is strongly connected to the amount of cold gas 
supply. This supports Mgas as the key ingredient to explain the mu¬ 
tual evolution of star formation and AGN accretion activity (e.g. 
Santini et al. 2012; Mullaney et al. 2012a; Rosario et al. 2013; 
local AGN review by Heckman & Best 2014). However, the link 
between BHAR and SFR that we observe, suggests that AGN ac¬ 
cretion is more closely related to the amount of cold star-forming 
gas (Mgas.sp), rather than to the total Mgas. Though in this work 
we are not able to distinguish between nuclear and global star for¬ 
mation, we note that if AGN accretion were tracing global cold gas 
mass, the known strong differences in star-forming efficiency (or in 
“depletion time”) between normal galaxies and luminous starbursts 
(e.g. Solomon et al. 1997; Saintonge et al. 2012) would predict a 
decreasing trend in BHAR/SFR with MS offset that we do not ob¬ 
serve (section 5.2.3). 

As argued in Section 5.2.2, our analysis may support the SFR 
as the original driver of the correlation with BHAR, though not well 
discernible from M, within individual redshift slices. By combin¬ 
ing all redshift bins, we are able to increase the statistics and infer 
more general trends between AGN accretion and galaxy proper¬ 
ties. Fig. 7 shows the sub-sample of our Herschel sources with FIR 
flux larger than the 80 per cent completeness threshold (see Section 
5.2.3). The middle and bottom panels in Fig. 7 show that the over¬ 
all correlation of BHAR with SFR is much narrower and significant 
than with M,. 

We proceed to a comparison between the BHAR-SFR rela¬ 
tion and the volume-averaged cosmic star-formation history and 
SMBH accretion history. As discussed in Madau & Dickinson 
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(2014, see their Fig. 15), both the global SFR density (SFRD) and 
black hole accretion rate density (BFIAD) peak around z~2, de¬ 
clining towards the local Universe. Despite the differences seen 
in recent derivations of accretion histories (e.g. Shankar et al. 
2009; Aird et al. 2010; Delvecchio et al. 2014; Ueda et al. 2014), 
all of them systematically support a slightly faster decay of the 
BFIAD since z~2 down to z~0, compared to the cosmic evolu¬ 
tion of SFH. By fitting the BHAD/SFRD ratio in the log-log space 
at 0<z<2, assuming the SFFI from Madau & Dickinson (2014) 
and an average BHAD evolution from the afore-mentioned liter¬ 
ature, we obtain BHAD oc SFRD^ '^^®'^. This relation is consis¬ 
tent with the slope evolution that we found between black hole 
and SFRs for our IR-selected sample: BHAR oc SFR'^ ®^°'^ (see 
Fig. 7, middle panel), supporting a scenario where the SMBH 
growth since z~2 follows a faster fading compared to their host 
galaxies. This redshift evolution is parametrized as BHAR/SFR 
oc (1-l-z)'^ ®^^'^ that we found through a linear regression fit in 
the log-log space. By integrating our estimates of average BHAR 
and SFR (see middle panel of Fig. 7) over cosmic time yields to 
(^/;r2°BHAR(z') dz')/(^/;=°SFR(z') dz') « 1/(3000±1500), 
marginally consistent with the prediction from the local Mbh- 
Mbuige relation (around 1/1000, e.g. Magorrian et al. 1998). This 
slight difference has been found and extensively argumented in var¬ 
ious studies (e.g. Jahnke et al. 2009; Rafferty et al. 2011; M12; 
C13). One possible factor is that our SFRs for main-sequence 
galaxies refer to the entire galaxy, including a large disk contri¬ 
bution, while only stars ultimately ending in the bulge should be 
counted in comparisons to the local Mbh-Mbuige relation. Sec¬ 
ond, since in the local Universe the black hole accretion is less 
radiatively efficient compared to non-local AGN (Merloni & Heinz 
2008), and the most powerful AGN are found in massive and pas¬ 
sively evolving systems, a sample of radiatively efficient AGN 
might miss a not negligible fraction of the local BH accretion rate 
density, leading to a smaller BHAR/SFR ratio. Finally, heavily ob¬ 
scured AGN that are not detected by deep Chandra surveys (e.g. 
Donley et al. 2012) may contribute to the cosmic AGN accretion 
history. 

As described in Section 5.2.3, in none of our redshift slices 
did we find significant evidence for a correlation between BHAR 
and MS offset. While, at fixed M*, outliers above the MS show 
enhanced BHAR (Fig. 5), this is washed out in a larger sample 
such as ours, that includes MS and above-MS bins of similar SFR. 

Recent studies based on deep Hubble Space Telescope (HST) 
images, such as the Cosmic Assembly Near-infrared Deep Extra- 
galactic Legacy Survey (CANDELS), allow thorough morpholog¬ 
ical analyses of star-forming galaxies at various redshifts. No dif¬ 
ference in morphology between AGN hosts and “inactive” galaxies 
has been found at 0.7 <z<3 (Cisternas et al. 2011; Schawinski et al. 
2011; Kocevski et al. 2011, 2012; Villforth et al. 2014), which re¬ 
shapes the relevance of major mergers and supports a less violent 
picture, where secular processes (e.g clumpy and/or unstable disks, 
Bournaud et al. 2012) play a major role in triggering both AGN 
and star formation activity. Observations carried out with Herschel 
(e.g. Santini et al. 2012; Mullaney et al. 2012b; Rosario et al. 
2013) have shown that the EIR-based sSFRs of moderately lumi¬ 
nous (Lx <10'*"^ erg s”'^) X-ray selected AGN hosts up to z~3 are 
almost indistinguishable from those of inactive galaxies. According 
to these studies, it is plausible to find a non-trend between offset 
from the MS and AGN activity in our Herschel-selected sample. 
The absence of a clear correlation is probably due to the fact that 
most of the luminous (Lx >10'*'^ erg s“^) AGN hosts are missed 


in our study due to the relatively small comoving volumes covered 
by our fields at low redshift. 

A further question concerns how the BHAR/SER ratio relates 
to major mergers. As noted, both morphology and total star for¬ 
mation of AGN hosts do not show evidence for a dominant role 
of major galaxy mergers in AGN triggering, of the type assumed 
by classical merger scenarios for the triggering of (luminous) AGN 
(Sanders et al. 1988; Earrah et al. 2001; Springel et al. 2005; Hop¬ 
kins et al. 2006; Sijacki et al. 2007; Di Matteo et al. 2008). If out¬ 
liers above the main sequence are mostly mergers, then questions 
on BHAR and BHAR/SER of mergers link back to the discussion 
of BHAR at different positions with respect to the main sequence. 
Over our full M* and SER range, we have found that MS offset 
is not a good predictor of absolute BHAR. However, at fixed M*, 
the enhanced SERs of above MS galaxies go along with enhanced 
BHAR, consistent with a BHAR/SER ratio that does not change 
with position with respect to the main-sequence. While the merger 
enhances both SFR and BHAR, it does so (on average) at a ratio 
consistent with that of MS galaxies. 

In summary, our results suggest that accretion onto the black 
hole and star formation broadly and on average trace each other, ir¬ 
respective of whether the galaxy is evolving steadily on the MS or 
bursting. This picture supports a causal connection between (radia¬ 
tively efficient) AGN activity and amount of the cold star-forming 
gas. Because of the different spatial and variability scales of the two 
phenomena, this connection is apparent only in sample averages. 


7 CONCLUSIONS 

We analysed the average BHAR in the SFR-M* plane at 0<z<2.5 
and investigated for the first time the mutual relations between 
BHAR and SFR, M, and sSFR in about 8600 Herschel-sdlscXed 
galaxies taken from the GOODS and COSMOS fields. The large 
statistics, together with the wealth of multi-wavelength data avail¬ 
able in these fields, allow us to explore a wide variety of star¬ 
forming galaxies and to characterize their individual SEDs through 
broad-band SED-fitting decomposition. Average AGN bolometric 
luminosities and accretion rates have been derived from Chan¬ 
dra X-ray observations, for both X-ray detected and undetected 
sources, under reasonable assumptions and widely used scaling fac¬ 
tors. 

Our main conclusions are as follows. 

(i) The average BHAR inferred in Herschel-selected galaxies 
shows positive evolution as a function of both SER and M* at 
z<0.8, while at higher redshift our data establish with >2a signif¬ 
icance the SER as the best predictor of AGN accretion. This may 
suggest that the galaxy SER is the original driver of the correlation 
with accretion onto the central black hole. 

(ii) Given the relation between BHAR and SFR, we found that 
the BHAR/SFR ratio slightly evolves with redshift, in agreement 
with a faster decline of the cosmic black hole accretion history with 
respect to the star formation history from z~2 to recent epochs, 
as published in the recent literature (e.g. Fig. 15 of Madau & 
Dickinson 2014). This evolutionary trend also leads to a lower 
BHAR/SFR ratio, albeit with a large associated error, compared 
to the predictions of the local Mbh-Mbuige relation. 

(iii) We compared the observed correlations between AGN ac¬ 
cretion and key galaxy properties with those presented by M12 and 
C13. Our analysis suggests that the (Mbh)-M* correlation at z~2 
claimed by Ml2 is likely a consequence of the trend with SFR and 
of the main-sequence relation that holds between the two. 
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(iv) These evolutionary trends of BHAR are plausible in the 
context of current studies of the evolution of the cold gas content 
in galaxies, if BHAR is on average linked to the content of dense 
star-forming gas. 


ACKNOWLEDGMENTS 

The authors are grateful to the referee, J. R. Mullaney, for his con¬ 
structive report. 

This paper uses data from HerscheVs photometers PACS and 
SPIRE. PACS has been developed by a consortium of institutes 
led by MPE (Germany) and including: UVIE (Austria); KU Leu¬ 
ven, CSL, IMEC (Belgium); CEA, LAM (Erance); MPIA (Ger¬ 
many); INAF-IESI/OAA/OAP/OAT, LENS, SISSA (Italy) and lAC 
(Spain). This development has been supported by the funding agen¬ 
cies BMVIT (Austria), ESA-PRODEX (Belgium), CEA/CNES 
(Erance), DLR (Germany), ASI/INAF (Italy), and CICYT/MCYT 
(Spain). SPIRE has been developed by a consortium of insti¬ 
tutes led by Cardiff Univ. (UK) and including: Univ.Lethbridge 
(Canada); NAOC (China); CEA, LAM (France); IFSI, Univ. Padua 
(Italy); lAC (Spain); Stockholm Observatory (Sweden); Imperial 
College London, RAL, UCL-MSSL, UKATC, Univ. Sussex (UK); 
and Caltech, JPL, NHSC, Univ. Colorado (USA). This devel¬ 
opment has been supported by national funding agencies: CSA 
(Canada); NAOC (China); CEA, CNES, CNRS (France); ASI 
(Italy); MCINN (Spain); SNSB (Sweden); STFC, UKSA (UK); and 
NASA (USA). 

Funding for SDSS-III has been provided by the Alfred P. 
Sloan Foundation, the Participating Institutions, the National Sci¬ 
ence Foundation, and the U.S. Department of Energy Office of 
Science. The SDSS-III web site is http : / /www . sdss3 . org/. 
SDSS-III is managed by the Astrophysical Research Consortium 
for the Participating Institutions of the SDSS-III Collaboration 
including the University of Arizona, the Brazilian Participation 
Group, Brookhaven National Laboratory, Carnegie Mellon Univer¬ 
sity, University of Florida, the French Participation Group, the Ger¬ 
man Participation Group, Harvard University, the Institute de As- 
trofisica de Canarias, the Michigan State/Notre Dame/JINA Par¬ 
ticipation Group, Johns Hopkins University, Lawrence Berkeley 
National Laboratory, Max Planck Institute for Astrophysics, Max 
Planck Institute for Extraterrestrial Physics, New Mexico State 
University, New York University, Ohio State University, Pennsyl¬ 
vania State University, University of Portsmouth, Princeton Uni¬ 
versity, the Spanish Participation Group, University of Tokyo, Uni¬ 
versity of Utah, Vanderbilt University, University of Virginia, Uni¬ 
versity of Washington, and Yale University. 

Funding for PRIMUS is provided by NSF (AST-0607701, 
AST-0908246, AST-0908442, AST-0908354) and NASA (Spitzer- 
1356708, 08-ADP08-00I9, NNX09AC95G). 

This research has made use of software provided by the Chan¬ 
dra X-ray Center (CXC) in the application package CIAO (v. 4.6). 

ID is grateful to Jillian Scudder for useful comments, and 
to Simonetta Puccetti and Nico Cappelluti for providing with X- 
ray sensitivity maps in the COSMOS field. ID is also thankful 
to Takamitsu Miyaji for his kind support in using CSTACK, and 
warmly thanks Fabio Vito for helpful suggestions with X-ray tools 
and routines to compute intrinsic X-ray luminosities. MB and 
GL acknowledge support from the FP7 Career Integration Grant 
“eEASy” (CIG32I9I3). 


References 

Ahn C. P. et al., 2014, ApJS, 211, 17 
Aird J. et al., 2012, ApJ, 746, 90 
Aird J. et al., 2010, MNRAS, 401, 2531 
Alexander D. M. et al., 2003, AJ, 126, 539 
Alexander D. M., Hickox R. C., 2012, New A Rev., 56, 93 
Amaud K. A., 1996, in Astronomical Society of the Pacific Con¬ 
ference Series, Vol. 101, Astronomical Data Analysis Software 
and Systems V, Jacoby G. H., Barnes J., eds., p. 17 
Barger A. J., Cowie L. L., Wang W.-H., 2008, ApJ, 689, 687 
Berta S. et al., 2013, A&A, 551, AlOO 
Berta S. et al., 2010, A&A, 518, L30 
Berta S. et al., 2011, A&A, 532, A49 
Bongiorno A. et al., 2012, MNRAS, 427, 3103 
Bournaud R, Dekel A., Teyssier R., Cacciato M., Daddi E., Juneau 
S., Shankar R, 2011, ApJ, 741, L33 
Bournaud R. et al., 2012, ApJ, 757, 81 
Boyle B. J., Terlevich R. J., 1998, MNRAS, 293, L49 
Brammer G. B., van Dokkum P. G., Coppi P, 2008, ApJ, 686, 
1503 

Capak P. et al., 2007, ApJS, 172, 99 
Chabrier G., 2003, ApJ, 586, L133 
Chen C.-T. J. et al., 2013, ApJ, 773, 3 

Ciliegi R, Gruppioni C., McMahon R., Rowan-Robinson M., 
2001, Ap&SS, 276, 957 

Ciotti L., Ostriker J. R, Proga D., 2010, ApJ, 717, 708 
Cistemas M. et al., 2011, ApJ, 726, 57 
Civano R et al., 2012, ApJS, 201, 30 
Coil A. L. et al., 2011, ApJ, 741, 8 

da Cunha E., Chariot S., Elbaz D., 2008, MNRAS, 388, 1595 
Daddi E. et al., 2010a, ApJ, 713, 686 
Daddi E. et al., 2007, ApJ, 670, 156 
Daddi E. et al., 2010b, ApJ, 714, LI 18 
Dekel A. et al., 2009, Nature, 457, 451 
Delvecchio I. et al., 2014, MNRAS, 439, 2736 
Di Matteo T., Colberg J., Springel V., Hernquist L., Sijacki D., 
2008, ApJ, 676, 33 

Donley J. L. et al., 2012, ApJ, 748, 142 
Elbaz D. et al., 2011, A&A, 533, Al 19 

Elmegreen D. M., Elmegreen B. G., Ravindranath S., Coe D. A., 
2007, ApJ, 658, 763 
Elvis M. et al., 2009, ApJS, 184, 158 
Farrah D. et al., 2001, MNRAS, 326, 1333 
Farrah D. et al., 2012, ApJ, 745, 178 

Feltre A., Hatziminaoglou E., Fritz J., Franceschini A., 2012, MN¬ 
RAS, 426, 120 

Ferrarese L., 2002, ApJ, 578, 90 
Forster Schreiber N. M. et al., 2009, ApJ, 706, 1364 
Fritz J., Franceschini A., Hatziminaoglou E., 2006, MNRAS, 366, 
767 

Gao Y, Solomon P. M., 2004, ApJ, 606, 271 
Garcla-Burillo S., Combes F., Schinnerer E., Boone R, Hunt 
L. K., 2005, A&A, 441, 1011 

Geach J. E., Smail I., Moran S. M., MacArthur L. A., Lagos 
C. d. R, Edge A. C., 2011, ApJ, 730, L19 
Gebhardt K. et al., 2000, ApJ, 539, L13 
Genzel R. et al., 2010, MNRAS, 407, 2091 
Gilli R., Comastri A., Hasinger G., 2007, A&A, 463, 79 
Grazian A. et al., 2006, A&A, 449, 951 
Griffin M. J. et al., 2010, A&A, 518, L3 


Mapping AGN accretion in the SFR-M^ plane 15 


Gultekin K., Cackett E. M., Miller J. M., Di Matteo T., Markoff 
S., Richstone D. O., 2009, ApJ, 706, 404 
Harrison C. M. et al., 2012, ApJ, 760, L15 
Heckman T. M., Best P. N., 2014, ARA&A, 52, 589 
Hickox R. C. et al., 2009, ApJ, 696, 891 

Hickox R. C., Mullaney J. R., Alexander D. M., Chen C.-T. J., 
Civano F. M., Goulding A. D., Hainline K. N., 2014, ApJ, 782, 
9 

Hopkins P. R, Hemquist L., 2009, ApJ, 694, 599 
Hopkins P. F., Somerville R. S., Hemquist L., Cox T. J., Robertson 
B., Li Y., 2006, ApJ, 652, 864 
Ilbert O. et al., 2009, ApJ, 690, 1236 
Ilbert O. et al., 2010, ApJ, 709, 644 
Jahnke K. et al., 2009, ApJ, 706, L215 
Kennicutt, Jr. R. C., 1998, ApJ, 498, 541 
Kocevski D. D. et al., 2012, ApJ, 744, 148 
Kocevski D. D. et al., 2011, ApJ, 736, 38 
Lanzuisi G. et al., 2013, MNRAS, 431, 978 
Le Floc’h E. et al., 2009, ApJ, 703, 222 

Leroy A. K., Walter F, Brinks E., Bigiel F, de Blok W. J. G., 
Madore B., Thomley M. D., 2008, AJ, 136, 2782 
Lilly S. J. et al., 2009, ApJS, 184, 218 
Lilly S. J. et al., 2007, ApJS, 172, 70 
Lutz D. et al., 2010, ApJ, 712, 1287 
Lutz D. et al., 2011, A&A, 532, A90 
Lutz D. et al., 2008, ApJ, 684, 853 
Madau P, Dickinson M., 2014, ARA&A, 52, 415 
Magdis G. E. et al., 2013, A&A, 558, A136 
Magdis G. E., Rigopoulou D., Huang J.-S., Fazio G. G., 2010, 
MNRAS, 401, 1521 

Magnelli B., Elbaz D., Chary R. R., Dickinson M., Le Borgne D., 
Frayer D. T., Willmer C. N. A., 2011, A&A, 528, A35 
Magnelli B. et al., 2014, A&A, 561, A86 
Magnelli B. et al., 2013, A&A, 553, A132 
Magorrian J. et al., 1998, AJ, 115, 2285 
Mainieri V. et al., 2007, ApJS, 172, 368 

Marconi A., Risaliti G., Gilli R., Hunt L. K., Maiolino R., Salvati 
M., 2004, MNRAS, 351, 169 
Merloni A., Heinz S., 2008, MNRAS, 388, 1011 
Mineo S., Gilfanov M., Lehmer B. D., Morrison G. E., Sunyaev 
R., 2014, MNRAS, 437, 1698 

Mineo S., Gilfanov M., Sunyaev R., 2012a, MNRAS, 419, 2095 

Mineo S., Gilfanov M., Sunyaev R., 2012b, MNRAS, 426, 1870 

Mullaney J. R. et al., 2012a, ApJ, 753, L30 

Mullaney J. R. et al., 2012b, MNRAS, 419, 95 

Netzer H., 2009, MNRAS, 399, 1907 

Noeske K. G. et al., 2007, ApJ, 660, L43 

Nordon R. et al., 2012, ApJ, 745, 182 

Nordon R. et al., 2010, A&A, 518, L24 

Oliver S. J. et al., 2012, MNRAS, 424, 1614 

Page M. J. et al., 2012, Nature, 485, 213 

Piconcelli E., Jimenez-Bailon E., Guainazzi M., Schartel N., 
Rodrlguez-Pascual P. M., Santos-Lleo M., 2005, A&A, 432, 15 
Pilbratt G. L. et al., 2010, A&A, 518, LI 
Poglitsch A. et al., 2010, A&A, 518, L2 

Rafferty D. A., Brandt W. N., Alexander D. M., Xue Y. Q., Bauer 
F. E., Lehmer B. D., Luo B., Papovich C., 2011, ApJ, 742, 3 
Ranalli P, Comastri A., Setti G., 2003, A&A, 399, 39 
Rodighiero G. et al., 2011, ApJ, 739, L40 
Rosario D. J. et al., 2013, ApJ, 763, 59 
Rosario D. J. et al., 2012, A&A, 545, A45 
Roseboom 1. G. et al., 2012, MNRAS, 419, 2758 


Saintonge A. et al., 201 la, MNRAS, 415, 32 
Saintonge A. et al., 201 lb, MNRAS, 415,61 
Saintonge A. et al., 2012, ApJ, 758, 73 
Salvato M. et al., 2011, ApJ, 742, 61 

Sanders D. B., Soifer B. T., Elias J. H., Madore B. F, Matthews 
K., Neugebauer G., Scoville N. Z., 1988, ApJ, 325, 74 
Santini P. et al., 2009, A&A, 504, 751 
Santini P. et al., 2014, A&A, 562, A30 
Santini P. et al., 2012, A&A, 540, A109 
Santos-Sanz P. et al., 2012, A&A, 541, A92 
Schawinski K., Treister E., Urry C. M., Cardamone C. N., Sim¬ 
mons B., Yi S. K., 2011, ApJ, 727, L31 
Schmidt M., 1959, ApJ, 129, 243 
Schreiber C. et al., 2014, arxiv: 1409.5433 
Scoville N. et al., 2007, ApJS, 172, 1 

Shankar F, Weinberg D. H., Miralda-Escude J., 2009, ApJ, 690, 
20 

Shao L. et al., 2010, A&A, 518, L26 

Sijacki D., Springel V., Di Matteo T., Hemquist L., 2007, MN¬ 
RAS, 380, 877 

Solomon P. M., Downes D., Radford S. J. E., Barrett J. W., 1997, 
ApJ, 478, 144 

Speagle J. S., Steinhardt C. L., Capak P. L., Silverman J. D., 2014, 
ApJS, 214, 15 

Springel V., Di Matteo T., Hemquist L., 2005, MNRAS, 361, 776 

Sutherland W., Saunders W., 1992, MNRAS, 259, 413 

Symeonidis M. et al., 2014, MNRAS, 443, 3728 

Tacconi L. J. et al., 2010, Nature, 463, 781 

Tacconi L. J. et al., 2013, ApJ, 768, 74 

Tozzi P. et al., 2006, A&A, 451, 457 

Trump J. R. et al., 2009, ApJ, 696, 1195 

Ueda Y, Akiyama M., Hasinger G., Miyaji T., Watson M. G., 
2014, ApJ, 786, 104 

Vasudevan R. V., Fabian A. C., 2007, MNRAS, 381, 1235 
Vattakunnel S. et al., 2012, MNRAS, 420, 2190 
Veilleux S. et al., 2009, ApJ, 701, 587 
Villforth C. et al., 2014, MNRAS, 439, 3342 
Vito F et al., 2014, MNRAS, 441, 1059 

Whitaker K. E., van Dokkum P. G., Brammer G., Franx M., 2012, 
ApJ, 754, L29 

Wuyts S. et al., 201 la, ApJ, 738, 106 
Wuyts S. et al., 201 lb, ApJ, 742, 96 
Xue Y. Q. et al., 2010, ApJ, 720, 368 
Xue Y. Q. et al., 2011, ApJS, 195, 10 


APPENDIX A: DETAILS ON X-RAY STACKING 

The stacking analysis allows to increase the signal-to-noise of in¬ 
dividually non-detected sources by grouping them together and 
piling-up their photon counts. This technique returns an average 
net (i.e. background subtracted) count rate for a given input list. 
As mentioned in Section 4.1, we stacked individual Herschel- 
selected galaxies not detected in X-rays. We masked all point-like 
sources listed in public X-ray catalogues (Alexander et al. 2003 for 
GOODS-N, Xue et al. 2011 for GOODS-S and Elvis et al. 2009 for 
C-COSMOS). We have assumed a circular shape for each mask, 
whose radius is 50 per cent larger then the area enclosing 90 per 
cent of photon counts from the source. The size of the Point Spread 
Function (PSF) follows a radial (and energy-dependent) profile in 
the Chandra Deep Fields, while in COSMOS it shows fluctuations 
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in between overlapping pointings (see Elvis et al. 2009). We ac¬ 
counted for that and ensured that the mean background level, after 
masking all X-ray sources, was fully consistent in each field with 
that presented in the corresponding reference papers. 

To optimize the final signal-to-noise ratio (SNR), we have 
taken a fixed source aperture (i.e. 2 arcsec radius) centered on each 
optical position. The photon counts have been corrected for off-axis 
angle and observed frequency^®. 

Each background region is an annulus centered on the cor¬ 
responding optical position, with inner radius placed 5 pixels (1 
pixel = 0.492 arcsec) apart from the PSF size enclosing 90 per cent 
source counts. The outer radius is 5 pixels apart from the inner one. 
This avoids contamination from the PSF of the stacked source. 

We subtracted the background counts, rescaled to the source 
area, from the photon counts collected within the source region, to 
get the net source counts. The exposure time is taken from energy- 
dependent time-maps, which return for each (x,y) position in the 
sky the effective exposure time (i.e. corrected for vignetting, bad 
pixels, dithering and including also a spatially dependent quantum 
efficiency). Finally, the average count rate is computed by summing 
(over all stacked sources) all net source counts and dividing by the 
total effective exposure time. 


APPENDIX B: STATISTICS IN THE SER-M* PLANE 

In Table B1 we show the values of BHAR and its uncertainties, as 
a function of SFR, M* and redshift, in order to quantify the colour- 
coded scale shown in Fig. 3. In addition, we list the number of 
sources included in each bin (total numbers. X-ray detected and 
X-ray undetected). Note that bins of SFR change as a function of 
redshift, as they are arranged to follow the cosmic evolution of the 
main-sequence relation. 


This con'ection has been made using the psf module (http : / /cxc . 
harvard . edu/ciao/ahelp/psf . html) from the Chandra Interac¬ 
tive Analysis of Observations (CIAO) package. 
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Table Bl. List of average BHAR, as a function of SFR, M* and redshift, with related Icr uncertainties. Alongside of each redshift range, we report the median 
redshift value ( 2 ) and the sSFRms that follows the main-sequence at ( 2 ). In some bins only an upper limit (upp) or a lower limit (low) is available. Numbers 
between parentheses in each bin represent: total number of Herschel sources, X-ray detected and X-ray stacked, respectively. 
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